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^ ■ Abstract 



We have shown that the lattice-gas model with four repulsive pair interaction 



constants (corresponding to the four nearest coordination shells) on a simple 
q | cubic lattice is a sufficient model to describe the main types of cation ordering 

in relaxors. The phase diagram, obtained by the cluster variation method, 

X 

shows the sequences of transitions between the phases 1 : 2 — > 1 : 1 
disordered and 1 : 2 — > disordered. 



I. INTRODUCTION 



Ferroelectric and piezoelectric properties of so-called relaxors, i. e. some complex mixed- 
metal perovskites with general chemical formulae A(B' X B" i_ x )0 3 , depend on the state of 
order in the (B', B") sublattice It is well-known that such relaxor materials as 

Pb(B'i/2Nbi/2)03 (B'=In and Yb) are ferroelectric when disordered and antiferroelectric 
when ordered [f|. The correlation between the order- disorder tendency of (B', B") ions and 
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the electromechanical (e/m) coupling factor, an important parameter for fabrication of high- 
quality piezoelectric materials, is also known [[J. Thus the compositional ordering in the (B', 
B") sublattice might be essential for the different applications in which relaxor materials 
are nowadays used: the A=Pb-based relaxors are widely used in piezoelectric transducers 
for mostly medical diagnostic application, in actuators, capacitors, ultrasound and sonar lis- 
tening devices, and the A=Ba-based relaxors are used in production of dielectric resonators 
and other high-frequency applications. 

Apart from well-known examples, such as BaTi03 with occupation of the (B', B") sub- 
lattice exclusively by 4-valent Ti ions (x = 0), there are two other main stoichiometrics of 
this sublattice with x = 1/2 and x = 1/3. While tetravalent (B', B") ions show no ordering 
(e. g. Pb(Zr a ,Ti 1 _ x )0 3 0), samples made of heterovalent ions, i.e. 3-valent B'= In, Sc, Yb 
and 5-valent B"= Nb, Ta or 2-valent B'= Zn, Mg, Ni, Co and 6-valent B"= W, Mo, in most 
cases exhibit long range order below an ordering transition temperature for x = 1/2. The 
same is true for most of the 2-valent B' and 5-valent B" ions at x — 1/3 (for experimental 
data on ordering types and temperatures see Tables in Refs [§,0). 

Basically, two types of (B', B") sublattice ordering are found experimentally. For x = 1/2 
compounds the 1:1 stacking of B' and B" layers perpendicular to the [111] axis of the simple 
cubic cell is observed in all A=Ba- and Pb-based A(B' x B"i_ x )03 materials. For x = 1/3 
the 1:2 structure with stacking of B', B", B" layers is found in all A=Ba-based relaxors 
with B'=Zn, Mg and B"=Nb, Ta. This type of ordering is the ground state for x = 1/3 
in a Coulomb-type of B'-B" interaction model 0. However the same four materials with 
Ba fully substituted by Pb in the A sublattice demonstrate either short-range 1:1 order- 
ing (e.g. Pb(Zn 1/3 Nb2/3)0 3 (PZN), Pb(Mg 1/3 Nb2/3)0 3 (PMN)) or long-range 1:1 ordering 
(Pb(Mg!/ 3 Ta2/ 3 )0 3 (PMT)). This structure might be viewed as an alteration of B', B" se- 
quence along the [111] axis with the B" site occupied by Ta(Nb) and B' by a disordered mix- 
ture of 2/3Mg(Zn)+l/3Nb(Ta) ||. Especially interesting is the type of (B', B") ordering ob- 
served in so-called relaxor ferroelectric alloys, for example (1— y)Ba(Zn 1 / 3 Ta2/ 3 )0 3 +?/BaZr0 3 
(BZT-BZ) or (1 - y)Ba(Mg 1/3 Nb 2 / 3 )0 3 +yBaZr0 3 (BMN-BZ). For small values of y < 0.02 
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the 1:2 ordering is observed. At intermediate values (0.04 < y < 0.25 for BZT-BZ and 
0.05 < y < 0.15 for BMN-BZ) the 1:1 ordering is observed in small microregions in a 
disordered matrix. At even higher y disordered structures are observed. 

Therefore it is interesting and important to find the appropriate model and from this 
model obtain the phase diagram (PD) which would account for known experimental facts 
including the disordered phase at high temperatures and the 1:1 and 1:2 phases at their 
corresponding stoichiometries. This problem is closely related to the problem of finding a 
minimal set of interaction constants in an Ising-type model for a simple cubic lattice. 

One of the first models, proposed to simulate the images found in HRTEM experiments 
for PMN, accounted for the nearest neighbor (NN) and next-nearest neighbor (NNN) in- 
teractions between (B', B") ions supplemented by an electrostatic (Coulomb) energy term 
E e i ||. At E e i = this Monte Carlo simulation yields the separation of initial 1:2 system into 
locally ordered negative 1:1 (Pb(Mgi/ 2 Nbi/2)03) 1_ domains surrounded by positive regions 
of (PbNbOa) 1 " 1 ". The minimization of electrostatic energy between these domains, realized by 
the inclusion of E e i, leads to much smaller 1:1 regions percolating the lattice. The obtained 
texture, which gives the superlattice reflections much closer to those seen in experiments, 
might be interpreted as a vanishing of the two-phase segregation. Another Monte Carlo 
simulation f2|, used the Coulomb interaction as the only type of interaction energy present 
in an ideal system of ionic charges. The comparison of the results of this model to different 
types of long-range structures in relaxors with x = 1/2 and x = 1/3 supported very much 
the viewpoint that the Coulomb interaction between ions might be the driving force of the 
ordering. 

There are also known several studies of the ground states of a simple cubic lattice using 
the Ising model with just few short-range interaction constants taken into account |9|-|TT||. 
However, the range of interaction of these models was not sufficient for the occurence of 
the ground states found in relaxors, the 1:2 phase in particular. Namely, the number of 
interaction constants were limited either to those three corresponding to the three main 
distances of the cube (i. e. NN, NNN and third- nearest neighbor (3NN)) PJlOl or to a very 



particular set of NN, NNN and fourth-nearest neighbour (4NN) interactions as well as to the 



smallest triangle interaction and linear triplet, but omitting 3NN interaction |TT|. Therefore 
much more important to the studies of stability of ordered phases in relaxors were the works 
of Burton [|3] and McCormack and Burton |12|] where the stability conditions for the 1:2 
phase were investigated in particular. These authors have concluded that using simple pair 
interactions in the Ising-type Hamiltonian up to fourth-nearest neighbors is not sufficient 
to obtain the experimentally observed superstructures as ground states. They obtained the 
PD with phase transition sequences 1 : 2 — > 1 : 1 — ■> disordered and 1 : 2 — > disordered 
using many-body interactions such as linear triplet and, so-called, cube-222. However the 
values of these interaction constants, chosen |J for the calculation of the PD, cannot be 
reconciled with the relative energies for several of the most probable phases obtained from 



first-principles calculation by the same authors [T2| . 

In this work we present the PD of a lattice-gas model with four repulsive pair interaction 
constants (NN, NNN, 3NN and 4NN) taken into account in a simple cubic lattice. Our calcu- 
lations, performed by the cluster variation method (CVM), show that this model is sufficient 
to explain the main ordering features of relaxors, and that no multi-point interactions are 
really necessary. 

II. MODEL AND PHASE DIAGRAMS 

We consider a lattice-gas model which is described by the Hamiltonian 

where rii = 1 when the site is occupied by the atom B' and when it is occupied by B". 



It can be shown ||T3| that the effective interaction constant Vij = v B ' B ' + v B " B " — 2v B ' B " , 



where vf? is the interaction energy between the atoms M and N, and /i = — fiB" + 
v b"b" _ v b'b" (yMN _ ^ vf- N ) is the chemical potential of the system expressed via the 
chemical potentials of the different cations. The lattice gas model can be mapped to the 
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more common Ising model in an external field via the relation between n, and the Ising spin 
variables Sj = ±1, = 2rij — 1, which implies the simple relation Vij = 4Jij between the 
Ising interaction constants Jjj and v^j. 

We are going to show that repulsive pair interactions in the four nearest coordination 
shells are sufficient to describe the experimentally observed 1:1 and 1:2 phases and the 
phase transitions between these phases. This statement contradicts earlier ground-state 
analysis jl2| , thus the Appendix is devoted to show that our approach is justified. Here we 



just want to give a simple indirect argument supporting our choice. 

According to the ground-state analysis [l(J for repulsive Vi, V2, t>3 > 0, v 4 = and the case 
^3 < 1*2/2, v 2 < t>i/4 + v 3 , four ordered ground-states take place at x = (n,) = 1/2, 3/8, 1/4 
and 1/8. The ground-state structure is 5(0, 6, 0; |) or 1:1 at x = 1/2, 5(0, 4, 0, |) (hereafter 
labelled as 3/8) at x = 3/8, 5(0,0,4 ±) (1/4) at x = 1/4 and 5(0, 0, 0; ±) (1/8) at x = 



1/8. The structure notations here, S(pi,p2,P3',x), are taken from Ref |T(J, and pi means 
the number of B'-B' bonds per B' atom in the z-th coordination shell of the B' X .B" 1 _ X 
structure (i=l, 2, 3 corresponds to NN, NNN and 3NN respectively). The phases are shown 
schematically in Fig .[I]. The 1/4 phase is bcc with a size of a new unit cell a = 2ao, the 1/8 
phase is sc with a = 2ao, and the 3/8 phase could be regarded as 1:1 phase with substracted 
1/8 structure (atoms are removed from the sites which are part of a = 2clq cubic structure). 
Here a is the lattice parameter of the initial sc unit cell. The ground-state at x — 1/3 
(the stoichiometric value for the 1:2 structure) is a mixture of 1/4 and 3/8 structures taken 
with proper weights. One can easily check that for v 4 = the energy of the 1:2 structure 
at x = 1/3 is exactly the same as that of the mixture, i.e. the 1:2 phase is just marginally 
unstable. By introducing an arbitrarily small v 4 > one can make the 1:2 structure stable 
albeit only at very low temperatures. This is because the energy of the mixture increases 
with f 4 , but the energy of the 1:2 phase is independent of v±. 

When choosing the empiric values of the interaction constants we assume V\ — 1, since 
the negative value would cause decomposition rather than ordering and the absolute value 
just gives the energy scale which can always be adjusted to the experiment. It would be 



interesting to explore all possible situations in the remaining three-dimensional space of v 2, V3 
and t>4 for which the ground state is the 1:2 structure at x — 1/3. This is too ambitious 
however, since not much is known about the PDs for the simple cubic lattice even for V4 = 0. 
Thus to simplify our task, we base our choice on the results of Ref. 0, that the interactions 
might be of electrostatic nature (similar atoms tend to repel each other) and take all Vi to 
be positive. 



The calculations of the PDs have been performed using the CVM ||14|| , a powerful vari- 
ational method to treat phase transition problems. In this method the local interactions 
between a particular group of spins (cluster) are exactly taken into account. The interaction 
of the cluster spins with those being outside the cluster are accounted for by the, so-called, 
effective fields which are the variational parameters in the version [15] of the CVM used 



here. The PD for t> 2 = 0.2, v 3 = 0.06 and v 4 = obtained using the 8+7 - point (cube + 
octahedron) approximation of the CVM is presented in Fig. |2|. The values of interactions 
were chosen to obtain the 1:1 phase down to x < 1/3 at intermediate temperatures (which 
are below the maximum of the 1:1 phase at x = 1/3, but above all other ordered phases). 
This PD could describe those A(B' 1 / 3 B" 2 /3)03 materials, which do not exhibit 1:2 ordering, 
but the ground state in this case would be a mixture of 1/4 and 3/8 phases. By setting 
V2 = V3 = we could make other structures disappear arriving at the classical antiferromag- 
netic Ising result with only 1:1 ordering. As was already mentioned, there is no 1:2 phase in 
this PD, since at x = 1/3 the energy of the 1:2 structure is equal to that of the proportional 
mixture of 1/4 and 3/8 phases. 

For v 4 > the 1:2 phase becomes energetically favourable with respect to the mixture 
of neighboring phases. The higher V4, the higher temperatures the 1:2 phase would reach. 
Other consequences of v 4 > on the PD are as follows. First, the domain of the 1:1 phase 
at intermediate temperatures should shrink at higher x, since the inclusion of V4 increases 
the energy of the 1:1 phase. Second, the degeneracy of the 1/8 and 3/8 phases is removed 
at v 4 > (at v 4 = any [100] or equivalent plane of the 1/8 phase can be moved along any 
vector in that plane and the depleted [100] planes of the 3/8 phase can be moved along i> 2 in 
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that plane without change of energy). Actually, for small v 4 > the ground state at x = 3/8 
is a new 3/8 phase (3/8'), which is obtained from the "original" 3/8 phase by shifting every 
second depleted [100] plane along the (Oil) direction (see Fig.^). Also at x — 1/8 the ground 
state for v± > is different from the "original" 1/8 phase (we do not examine this phase 
as it has no relevance to the vicinity of x = 1/3). Third, a number of ordered structures 
with longer period should appear below x = 1/8 and above x = 3/8, which we also do not 
consider here. For small V4 > 0, the transition from 1:2 to 3/8' phase takes place between 
x = 1/3 and x = 3/8. Another transition from 1/4 to 1:2 occurs between x = 1/4 and 
x = 1/3. The PD for v-i = 0.2,^3 = 0.06 and V4 = 0.02 is presented in FigfJ. The PD is 
not complete as the phases to the left of the 1/4 phase and to the right of the 3/8' phase 
are not considered. We obtain a rather narrow 1:2 phase in (T,x) coordinates. Meanwhile 
in the (T, //) diagram, there is quite a large domain of chemical potential corresponding to 
the 1:2 structure. It should be also noted that for V4 = 0.02 we could not obtain the 3/8' 
phase by the CVM, since the mixture of the 1:2 and off-stoichiometric 1:1 phase has lower 
energy. We calculated the high-temperature metastability limit of the 3/8' phase by Monte 
Carlo simulations only (see Fig.[3]b). Still we can conclude that if the 3/8' phase exists, it 
exists only below T ps 0.08. 

For 2.01 < fj, < 2.06 and 0.13 < T < 0.25 the 8+7 - point approximation of the 
CVM has unexpectedly shown some instability resulting in the gap on the phase boundary 
between 1:1 and 1:2 phases and the strange reentrance of the second order transition line 
between disordered and 1:1 phase. It is interesting to note, that reducing the variational 
space of the CVM by omitting all the correlations with more than two sites, the problems 
related to the second order phase transition line disordered - 1:1 disappear (dotted line 
in Fig.|3|), but those related to the transition 1:1 - 1:2 become worse. Such a behavior 
of the 8+7 - point approximation is clearly related to the increased range of interactions 
which can not be properly accommodated within the cube and octahedron clusters. It 
seems that these problems appear at any v 4 > (we checked the values between 0.01 and 
0.04, also slightly varying V2 and V3). At v 4 = we had no similar problems with the 
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8+7 - point approximation. To be sure, we have performed calculations with a 27 - point 
(2ao x 2a x 2a cube) approximation of the CVM. Besides the increase of accuracy due to 
larger basic cluster (a general feature of the CVM), the 27-point cluster has an advantage 
of including all interactions into single basic cluster, thus avoiding occasional breakdowns, 
sometimes observed when basic clusters are too small. Unfortunately, the minimisation of 
the free energy then becomes more complicated, since the number of variational parameters 
drastically increases. Since it is not feasible to treat all 2 27 pa 134 mln configurations of the 
27 - point cluster, we excluded the possibility for the NN sites to be occupied {v\ = oo). 
This NN exclusion makes a big difference at high temperatures, as the 1:1 phase is always 
present above some critical concentration. However we expect only small changes at low 
temperatures, where 1:2 and the rest of ordered phases exist, since due to v\ ^> T the 
probability of the NN being occupied is very low, especially for small concentrations. The 
NN exclusion decreases the number of configurations of the 27 - point cluster down to 70663, 
which is tractable with some effort. The results of the 27 - point CVM are presented by 
dashed lines in Figs. ^ and |3|. Indeed, while there are large discrepancies between the results 
obtained by the 8+7- and 27 - point approximations at higher temperatures due to the 
exclusion, at low temperatures the lines are pretty close. It is rather difficult to compare 
accuracies of the two approximations, since formally they describe different models. The 
27 - point approximation does not suffer from the problems of the 8+7 - point one and 
confirms its results by filling the gap in the PD. To be sure, we have also applied the NN 
exclusion to the 8+7 - point approximation (to check whether the exclusion itself could hide 
the problems with the CVM), but the problems remained. 

To support our CVM calculations, we have performed some limited Monte Carlo simula- 
tions using a simple Metropolis algorithm on a 24 x 24 x 24 lattice with periodic boundary 
conditions and around 50000 sweeps for each point. These points are plotted in Fig.|3] with 
circles. First-order transitions are difficult to locate accurately by the Monte Carlo method, 
so the points on the first-order phase transition lines are actually metastability limits of the 
corresponding low-temperature phase, a true phase transition should be slightly below. It 
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is much easier to find metastability limits, if one starts from a perfectly ordered structure 
at low temperatures and gradually increases the temperature until the order changes or is 
destroyed. We show in Fig.^, that the metastability limits obtained by the CVM and Monte 
Carlo calculations are quite close, thus we expect the location of the true first-order phase 
transition to be correct too. It is seen that these results nicely confirms our CVM data. 

In summary, a lattice gas model for a simple cubic lattice with four nearest repulsive 
pair interaction constants was proposed. The presented PD posesses all features of cations 
ordering observed in relaxor compounds, namely, a phase transition sequence 1 : 2 — ^ 1 : 1 — > 
disordered at x — 1/3. The increase of v± from 0.02 to 0.03 shifts the phase boundary between 
disordered and 1:1 phases to the region x > 1/3, thus resulting at x = 1/3 in a phase 
transition 1 : 2 — ^ disordered, which was also observed in a number of relaxor materials. 
Another interesting feature of the PD is that a pure 1:2 phase reaches its stoichiometric 
concentration only at T = 0, while at any finite temperature it separates into a mixture of 
the 1:2 phase with a slightly lower concentration and a small amount of adjacent phase (e.g. 
1:1). This indicates that a slightest increase in concentration above stoichiometric value 
tends to destroy the 1:2 ordering, making it favorable for the excess ions to form islands 
of another higher concentration phases. Could this be an explanation of the nanoregions, 
which are observed in some A(B'i/ 3 B" 2/3)03 compounds and, most likely, are responsible 
for the relaxor properties, besides the electrostatic arguments? 

This work was supported by The Swedish Academy of Sciences, The Swedish Insti- 
tute, The Swedish Natural Research Council and SSF (Swedish Foundation for Strategic 
Research) . SL is indebted to The Royal Institute of Technology, Stockholm for kind hospi- 
tality. 



APPENDIX A: 



Here we apply the partial vertex enumeration (PVE) procedure 112] for interaction clus 



ters to our four two-point clusters corresponding to NN, NNN, 3NN and 4NN pairs. The 
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PVE technique for the ground state problem is based up on the more general vertex enu- 
meration of the configuration polyhedron method ||13| . It is easier described in terms of 
the Ising model and the outline is as follows. Consider a cluster a of \a\ lattice sites which 



includes all interactions of the Hamiltonian. It is shown |13| that the probability function 
pa{{ 8 }) for the configuration {s} (the set of Ising spins at every site i) of the cluster a is 
a linear function of the spin products sp = Yliep s i °f an subclusters j3 of the cluster 

f3Ca,P^0 

where up = (sp) are cluster correlation functions. The requirement ^ p ^ 1 imposes 
constraints for the allowed values of up, and equations 

Pads}) = (A2) 

for each symmetry distinct configuration {s} define hyperplanes in the space spanned by up. 
These hyperplanes are boundaries in a space of allowed values of vp which is a configurational 
polyhedron (CP). The possible ground states of the model reside in the vertices of the CP, 
because the energy is a linear function of up, and the extremum of the linear function is 
reached on the boundary of its domain of definition. Since not all up necessarily contribute 
to the energy E = (Tt) , one needs to project the CP on the subspace of those up, which 
contribute to the energy. 

However it may be impractical to make a full vertex enumeration for large clusters, when 



interactions beyond the 3NN are considered. Therefore the PVE is used (see Ref [|12]) to 
determine whether in the space of Up(f3 C a) included in a certain cluster a (or several 
clusters), the point, corresponding to the given structure, is a vertex of the CP. It is a vertex 
if the number of faces it belongs to is not less than the dimensionality of the space (number 
of distinct up). Then for a given (periodic) structure one obtains values of up by averaging 
over the unit cell of the structure, and checks how many equations ( |A2| ) are satisfied by 
these values. It is found (TJ| that at least some many-body interactions in the cube and a 



linear triplet are necessary for the 1:2 phase to be a vertex. 
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Consider PVE for our four pair clusters (a = 1,2,3,4). The configuration probabilities 
for each cluster are p a — |[1 + rn(si + Sj) + r a SiSj], or 

p a (±±) = ^(l±2m + r a ) ^0, p a (+-)=p a (-+) = ^(l-r a )^0. (A3) 

Here m = (sj) is a magnetization, r a = (siSj)i j ea is a pair-correlation function of the 
cluster a, and these five parameters span a 5-dimensional space. For the 1 : 2 structure 
m = r\ = r4 = —1/3, T2 = 1/3, r% = 0, and two configuration probabilities are zero, namely 
pi(++) and pi(++). This accounts for 2 faces of the 5-dimensional space. If there are no 
other constraints on the m,r a (a = 1, . . .,4), the structure cannot be a vertex. But there 
are other constraints. The point is, that p ^ are sufficient conditions only if all vp are 
defined on the same cluster. This is not the case here and neither it was in Ref fL2| . There 



are in general additional relations between the probabilities for configurations defined on 
different clusters, since one can imagine these clusters being parts of a single larger cluster. 
Let's take a simple example. Consider the 4NN pair (a = 4) and mark its sites with (+) 
if Si = +1, (— ) if Sj = — 1, and (•) for the site in the middle (to make it different from the 
NN pair). Also consider for a moment a linear triplet which contains a 4NN pair and two 
NN pairs. Since p(si,Sj) = ^2 Sk=±1 p(si, Sj, Sk) and all probabilities are non-negative, it is 

obvious that p(+ ■ — ) = p(+ H — ) + p(H ) ^ p(++) + p( )? which using ( [A3] ) gives a 

relation 

1 - r 4 ^ 2 + 2r x . (A4) 

This relation defines a new constraint and a new hyperplane, which is not equivalent to any 
of the constraints in ( |A3|) . It is not that simple to obtain all the constraints. In general one 
has to consider clusters large enough to contain all the interaction clusters (e.g. a double- 
cube), then to project the resulting CP on the space of m,r a (a = 1,...,4), and see if 
the vertex coresponding to the desired structure survives. Fortunately, here we can use the 
result [12] that the cube+4NN pair interaction cluster is missing one constraint for the 1:2 



phase to be a vertex. Since the constraint (|A4|) , which is not equivalent to those already 
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considered in Ref |12| , applies here also, and the corresponding hyperplane includes the 1:2 
phase, we conclude that the 1:2 structure is a possible ground state. 
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FIGURES 

FIG. 1. Schematic view of ordered phases mentioned in the text. The structures are projected 
onto the (001) plane. Open circles indicate a sequence of B" atoms along the [001] direction. 
Half-filled black-white (white-black) circles represent B', B" (B", B') atoms alternating along the 
[001] direction. 

FIG. 2. Phase diagrams (T,x) (a) and (T, //) (b) obtained by the CVM with the interaction 
constants Vi = 0.2, V3 = 0.06 and = 0. Solid line - 8+7 - point approximation (v\ = 1), dashed 
line - 27-point approximation [v\ = 00). The two-phase regions, marked by numbers are: 1 - 
disordered+1/4 and 2 - 1:1+3/8 

FIG. 3. Phase (T, x) (a) and (T,n) (b) diagrams obtained by the CVM with the interaction 
constants V2 = 0.2, V3 = 0.06, and V4 = 0.02. Solid and dashed lines - same as in Fig.||, dotted line 

- 8+7-point approximation with the variational space, reduced to one- and two-point correlation 
functions only. Monte Carlo results are shown by circles. Phases to the right of the 3/8' phase 
and to the left of the 1/4 phase are not shown. The two-phase regions, marked by numbers are: 1 

- disordered+1/4 and 2 - disordered+l:2 

FIG. 4. The metastability limits of 2:1 and 1:1 phases at fj, = 2.1 obtained with the CVM 
(solid lines) and Monte Carlo calculations (circles). The energy E = (TC) per site is presented for 
V2 = 0.2, V3 = 0.06, and V4 = 0.02. The dashed lines are CVM free-energies which intersect at the 
first-order phase transitions (pointed to by an arrow). 
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